Anomalies in water as obtained from computer simulations of the 
TIP4P/2005 model: density maxima, and density, isothermal 
compressibility and heat capacity minima 

Helena L. Pi, Juan L. Aragones, Carlos Vega*, Eva G. Noya, 
Jose L. F. Abascal, Miguel A. Gonzalez and Carl McBride 
Departamento de Qmmica Fisica, 
Facultad de Ciencias Qmmicas, 
Universidad Complutense de Madrid, 
28O4O Madrid, Spam 

(Dated: May 25, 2009) 

Abstract 

The so-called thermodynamic anomalies of water form an integral part of the peculiar behaviour 
of this both important and ubiquitous molecule. In this paper our aim is to establish whether the 
recently proposed TIP4P/2005 model is capable of reproducing a number of these anomalies. Using 
molecular dynamics simulations we investigate both the maximum in density and the minimum 
in the isothermal compressibility along a number of isobars. It is shown that the model correctly 
describes the decrease in the temperature of the density maximum with increasing pressure. At 
atmospheric pressure the model exhibits an additional minimum in density at a temperature of 
about 200K, in good agreement with recent experimental work on super-cooled confined water. 
The model also presents a minimum in the isothermal compressibility close to 310K. We have also 
investigated the atmospheric pressure isobar for three other water models; the SPC/E and TIP4P 
models also present a minimum in the isothermal compressibility, although at a considerably lower 
temperature than the experimental one. For the temperature range considered no such minimum 
is found for the TIP5P model. 
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I. INTRODUCTION 



The study of polar fluids has long been a topic of interest, even more-so since the devel- 
opment of perturbation theories and computer simulation techniquesii^. Amongst the many 
polar molecules, water stands out in particular. Water is a fascinating molecule, both from a 
practical and from a fundamental point of view. In the liquid phase water presents a number 
of anomalies when compared to other liquids^i^i^i^i^, whilst the solid phase exhibits a highly 
complex phase diagram, having at least flfteen different crystalline structures^i^. Due to 
its importance and its complexity, understanding the properties of water from a molecular 
point of view is of considerable interest and presents a veritable challenge. 

Computer simulations of water started with the pioneering papers by Watts and Barker- 
and by Rahman and Stillinger— about forty years ago. However, a key issue that still exists 
when performing simulations of water is the choice of model for the pair potential that is used 
to describe the interaction between moleculesii^^^^. The SFG^, SPC/Eii, TIP3P-i^, 
TIP4P^^ and TIP5P— models have become highly popular among the large community of 
people simulating water or water solutions. Each of the aforementioned models are rigid and 
non-polarisable, but naturally, real water is both flexible and can be polarised. It is almost 
needless to say that these models represent a very simplifled version of the true interactions 
between real water molecules. Water also exhibits important quantum nuclear effects, so a 
realistic description of water should also take this into account^^-. That said, it is of interest 
at this stage to analyse how far it is possible to go in describing real water using these simple 
models. 

In recent years the inexorable increase in computing power has permitted the calculation 
of properties that were previously inaccessible to simulations. These properties can be 
used as new 'target quantities' when fltting the parameters for any new potential. More 
importantly, some of these properties have provided stringent tests for the existing water 
models. In particular, some of the authors have recently determined the phase diagram 
for different water models and have found that the performance of the models can vary 
signiflcantly2i>22ii22i>2ii2^. For example, it has been found that TIP4P provides a qualitatively 
correct description of the global phase diagram of water, whereas the SPC, SPC/E, and 
TIP5P models do not.^^ For the SPC and SPC/E models the melting temperature and 
the maximum in density of water occur at temperatures far below the experimental values. 
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Taking this into account some of the authors proposed a new rigid non-polarisable model 
with the following constraints: that it should be based on the TIP4P model, since this model 
provides a reasonable description of the phase diagram of water. It should also reproduce 
the maximum in density of liquid water (notice that most of models fail when it comes to 
predicting the location of the density maximun>2i). Finally, the model should account for 
the vaporisation enthalpy of real water, but only after incorporating the self-polarisation 
energy correction proposed by Berendsen et al. for the SPC/E model^^. It was from these 
considerations that the TIP4P/2005 model of water arose.-^ In the paper in which the model 
was presented it was shown that TIP4P/2005 correctly describes the global phase diagram 
of water, the diffusion coefficient at atmospheric pressure and temperature, the maximum 
in density along the atmospheric pressure isobar, the density of the ice polymorphs, and the 
equation of state of liquid water at high pressures. Afterwards, additional studies have shown 
that the model is also able to provide a good description of the vapour-liquid equilibria^ 
and the surface tensiort^S. Further information concerning the performance of TIP4P/2005 
and its comparison with other water models can be found elsewhere^^. The overall results 
indicate that TIP4P/2005 is probably close to being the best rigid non-polarisable model 
that can be achieved; any significant improvement would require the introduction of 'new' 
physics, i.e., fiexibility, polarisability and nuclear quantum effects. 

In this paper we are interested in the thermodynamic response functions in the region 
where water exhibits 'anomalous' behaviour. In particular, these are the expansion coef- 
ficient (which vanishes at temperatures close to the melting point, resulting in the well 
known maximum of density at about 4°C at atmospheric pressure) and the isothermal 
compressibility^^ (^^^ ) which shows a minimum at 46.5°C at p=l bar. Moreover, we shall 
investigate whether the model is able to describe the pressure dependence of these ther- 
modynamic properties, namely, the decrease in the temperature of the density maximum 
and the shift towards slightly higher temperatures of the minimum in kt as the pressure 
increases. It will be shown that the model is indeed able to describe these two features 
quite well, which provides further evidence of the robust behaviour of the model even when 
estimating properties that were not taken into account in the original fitting process. 
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II. METHODOLOGY 



The interaction between water molecules will be described by the TIP4P/2005 model.— 
In this model, a Lennard- Jones centre is located on the oxygen atom, positive charges are 
placed on the hydrogen atoms and a negative charge is located at the site M situated along 
the H-O-H bisector. For the simulations we have used the molecular dynamics package 
GROMACS (version 3.3).'^^ The Lennard- Jones potential has been truncated at 9.0 A. Long 
range corrections were applied to the Lennard- Jones part of the potential (for both the 
energy and pressure). Ewald summations were used to deal with electrostatic contributions. 
The real part of the Coulombic potential was truncated at 9.0 A. The Fourier component 
of the Ewald sums was evaluated by using the Particle Mesh Ewald (PME) method of 
Essmann et al— The width of the mesh was 1 A and a fourth degree polynomial was used. 
The simulation box was cubic throughout the whole simulation and the geometry of the 
water molecules was enforced using constraints.— The temperature was set by using a 
Nose-Hoover^^'^^ thermostat with a relaxation time of 2 ps. To maintain constant pressure 
an isotropic Parrinello-Rahman barostat^"^ with a relaxation time of 2 ps was used. As a 
check, at two pressures, p=l bar and p=400 bar, Monte Carlo simulations were performed 
using a bespoke program. The Monte Carlo densities were in complete agreement with those 
obtained from molecular dynamics using GROMACS. 

To determine the maximum in density, molecular dynamics simulations have been per- 
formed along the isobars p=l, 400, 1000, and 1500 bar. The number of molecules used in 
the simulations was 256. Long runs are required to determine the maximum in density; for 
each thermodynamic state twenty million time-steps were performed. Since the time step 
was 2 fs, the results presented here are an average of the properties of the system obtained 
from runs of 40 ns. The simulations were started at high temperatures, and the final con- 
figuration of a particular run was used as the initial configuration for a lower temperature 
simulation. Typically about 6 to 8 different temperatures were studied along each isobar. 
The isothermal compressibility is defined as: 



The literature for the isothermal compressibility of water models is rather scant^i^. The 
computational overhead required for a study of the variation of kt with temperature with 




(1) 
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sufficient accuracy has, until now, been prohibitive. Here we have evaluated the isother- 
mal compressibility for the following water models: SPC/E^^, TIP4P— , TIP4P/2005^ and 
TIP5P— along the atmospheric pressure isobar. The isothermal compressibility has been 
evaluated using the well known fluctuation formula 

jV') - {Vf 
'^^ = kT{V) ■ 

The volume fluctuations were typically averaged over 20 ns using a time-step of 1 fs. Pre- 
liminary results for TIP4P/2005 model were sufficiently close to the experimental value as 
to warrant a more precise calculation. Thus, for this model, the isothermal compressibil- 
ity has been determined via two different procedures using, in both time-step of 
0.5 fs and a sample size of 500 molecules. The flrst method is the application of the volume 
fluctuation formula, typically over 40 million molecular dynamics time-steps, for a total 
simulation length of about 20ns. The relative uncertainty of the calculated compressibilities 
is about ±3%. The second route is the numerical evaluation of the derivative appearing 
in the deflnition of kt (Eq. [1]). To this end, the results for the equation of state of flve 
different state points (at different pressures but at the same temperature) were fltted to a 
second degree polynomial. In addition to the state point of interest, for which the volume is 
already known from the runs used to calculate kt via the fluctuation route, four additional 
runs are required: two runs at the same temperature but at pressures of 200 and 400 bar 
higher, and another two runs at pressures 200 and 400 bar lower. The simulation length of 
these additional runs was 6.5 ns , so in the end the computational cost of both routes is 
almost the same. We have also computed the heat capacity at constant pressure Cp for the 
TIP4P/2005 model. Cp is deflned as : 

Thus in order to compute Cp the enthalpy at each temperature was flrst calculated. These 
values for the enthalpy were then fltted to a polynomial and this flt was then differentiated 
with respect to temperature to obtain Cp{T). 



III. RESULTS 



The simulation results for the densities are reported in Tables [U HIl IIIII and HVl for p=l, 
400, 1000 and 1500 isobars, respectively. For the four isobars considered a clear maximum 
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in density has been found. The results for p=l bar and p =400 bar, along with experimental 
measurements, are plotted in Fig. [TJ It can be seen that the model is able to reproduce 
the experimental data quite nicely. The density maximum for each isobar was obtained by 
fitting the densities around the maximum (typically 5 or 6 densities were used in the fit) to 
a second or a third degree polynomial. The values of the maximum are reported in Table 
IVl In Fig. [2] the temperatures at which the maxima appear (TMD) are compared to the 
experimental values^. As can be seen, the agreement is rather good. The model correctly 
predicts a decrease of about 33 K in the temperature of the maximum when going from 
atmospheric pressure to a pressure of about 1500 bar. In other words, 45 bar are required 
to induce a decrease in the TMD of about one degree. For D2O Angell and Kanno found a 
similar lowering of the TMD with pressure.™ For the TIP4P/2005 not only does the TMD 
decrease with pressure but the melting temperature does too. In Figl2] the melting curve 
of ice I/i (taken from our previous work-2^'^'^) is plotted along the TMD curve determined 
in this work. As can be seen, both curves have a negative slope. Notice also that at room 
pressure the distance between the melting curve of ice I/, and the TMD is about 25K, which is 
in contrast with the 4K of difference found experimentally. The impossibility of reproducing 
simultaneously both the TMD and the melting temperature in models having three charges 
has been discussed previously2^>2I. 

Another interesting issue is the behaviour of the density along the atmospheric pressure 
isotherm at very low temperatures. Experimentally, it is not possible to obtain the density of 
water at temperatures below 233 K (the homogeneous stability limit of water at atmospheric 
pressure^) due to spontaneous nucleation and freezing. However, recently, it has been pos- 
sible to avoid the formation of ice by confining water in pores a few nanometers in diameter 
(most probably because of the decrease in the freezing point of water due to confinement as 
described qualitatively by the Gibbs-Thomson approach^). Thanks to this the density of 
deuterated water has recently been determined for the first time^ for temperatures as low 
as 160 K, reporting the existence of a minimum in the density at a temperature of about 
200 K. A similar study (for water instead of deuterated water) was performed by Mallamace 
et al^^. By using infra-red (IR) spectroscopy, it was possible to determine the density of 
liquid water in the super-cooled regime down to 150 K. In Fig. [3] the reported experimental 
densities (from IR spectroscopy of a sample within a 1.5 nm pore) for liquid water along 
the atmospheric pressure isobar are compared to those obtained in this work from molecular 



6 



dynamics simulations of the TIP4P/2005 water model. The agreement is surprisingly good, 
and the location of the minimum around 200 K is described properly by TIP4P/2005. Such 
a minimum in density has also been found^*^ in computer simulations of the TIP5P and 
ST2 models of water, but the agreement with experiment was not as good as that exhib- 
ited by T1P4P/2005. For temperatures below this minimum the density increases as the 
temperature decreases (as in a normal fluid). In Fig. H] we compare the equation of state 
of supercooled water to the equation of state of ice Ih for the TIP4P/2005 model along the 
room pressure isobar.— The minimum in density of supercooled water occurs just when the 
density is approaching that of the ice Ih- However, we did not succeed in obtaining ice Ih by 
cooling water; the radial distribution functions of the supercooled water (at room pressure) 
below 200K are clearly different from those of ice Ih- Rather the minimum in density corre- 
sponds to the formation of a glassy state. Note that the existence of such a minimum is not 
only restricted to water but is also present in materials such as telluriunt^. Interestingly 
the melting curve for tellurium exhibits re-entrance . Such re-entrant behaviour was also 
found in our studies of the phase diagram of water models^. 

We have also computed the self-diffusion coefficient at p =1 bar and p=1000 bar. The 
results, presented in Fig. [5], show that the diffusion coefficient drops significantly as the tem- 
perature decreases. The decay is less pronounced at higher pressures. At low temperatures 
the diffusion coefficient increases significantly with pressure (this behaviour is anomalous^^ 
since, for most fluids, it decreases with increasing pressure). One can imagine that the 
application of pressure somehow breaks the hydrogen bond network, thus aiding diffusion 
processes. Some indirect evidence of the decrease in hydrogen bonding with pressure can be 
obtained from the analysis of the different contributions (Lennard- Jones and Coulombic) to 
the residual internal energy. In Table |Vl] the different contributions to the residual internal 
energy obtained from molecular dynamics of the TIP4P/2005 model along the T=224.6K 
isotherm (at different pressures) are given. Notice the large decrease (in absolute value) of 
the Coulombic energy with pressure, clearly pointing to a reduction of hydrogen bonding 
with pressure. The Lennard- Jones contribution is more repulsive at low pressures (even 
though the density is lower) as a consequence of the stronger hydrogen bond found at lower 
pressures. At temperatures above 280 K the diffusion coefficients at the two pressures con- 
sidered become virtually identical as the differences fall to within the statistical uncertainty. 
The agreement with the experimental values of the diffusion coefficient^ is quite good. 
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Let us now focus on another of the 'anomalous' properties of water, the isothermal com- 
pressibility. The experimental measurements^i^i^i^ show that, at atmospheric pressure, 
kt drops as the temperature increases from the melting temperature up to 46.5°. Above this 
temperature water behaves as a normal liquid and the isothermal compressibility increases 
with temperature. It is also well established from experiment that the temperature for which 
hiT is minimal shifts slightly towards higher values as the pressure increases. Contrary to 
the maximum in density, the ability of the water models to account for the compressibility 
minimum has not yet been established. As mentioned in the previous section, simulation 
studies of the isothermal compressibility are few and far between^i^. Moreover, the com- 
putational resources available did not allow the extended simulations needed to calculate kt 
with the precision required to determine whether the most common water models predict 
the compressibility minimum. For this reason we have calculated the compressibility for a 
select few 'popular' water models. The results are presented in Table |VTT] for the SPC/E, 
TIP4P and TIP5P models and in Table Ivml for TIP4P/2005. The calculations for SPC/E, 
TIP4P and TIP5P were performed at atmospheric pressure for three different temperatures, 
namely 260 K, 300 K and 360 K. Fig. shows that the temperature dependence of the 
isothermal compressibility for TIP5P does not follow the experimental pattern since, for this 
model, kt is a monotonously increasing function along the whole experimental liquid range. 
Despite the fact that the slope of the TIP5P curve is opposite to the experimental one, the 
values of are coincident in a narrow range of temperatures because the simulation and the 
experimental curves cross at a temperature close to the freezing point of liquid water. For 
temperatures near the boiling point the TIP5P model fails completely, where the predicted 
compressibility at 360 K is almost twice the experimental value. 

As for the performance of the SPC/E and TIP4P models, both provide fairly similar 
results. In fact their curves are parallel, showing a more or less defined minimum around 
270 K (more computations would be needed to determine the precise location of the minima). 
The results of SPC/E are somewhat shifted to higher temperatures with respect to those 
of TIP4P, which results in a slightly better agreement with the experimental data. At high 
temperatures, the differences between simulation and experiment are quite noticeable for 
both SPC/E and TIP4P (though less dramatic than in the TIP5P case). 

As can be seen in Fig. |6]d, the TIP4P/2005 model provides an excellent description of 
the isothermal compressibility of water. The compressibilities obtained from the two routes 
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(differenciation of the equation of state and fluctuation formula) were found to be mutually 
consistent. The departures of the calculated compressibilities with respect to the experi- 
mental values^'^'^i^ are in general smaller than 7%. By fitting the TIP4P/2005 compress- 
ibilities at atmospheric pressure to a third degree polynomial, a minimum in the isothermal 
compressibility is found for a temperature of about 37°C which is in good agreement with 
the experimental value (46.5°C). Thus the model is able to accurately reproduce not only 
the isothermal compressibility along the whole liquid range at atmospheric pressure but also 
the compressibility minimum. This assertion is particularly true when one compares the 
TIP4P/2005 predictions with those of the other models, shown in Fig [6] (notice that we have 
plotted both panels using the same scale) . In Fig El we also present the compressibility 
results at a higher pressure (1000 bar). TIP4P/2005 predictions for this isobar are slightly 
better than those for atmospheric pressure. A compressibility minimum is also found at 1000 
bar, and, in accordance with the experiment, the minimum appears at an slightly higher 
temperature than it does at atmospheric pressure^. 

Finally we have examined the behaviour of the heat capacity at constant pressure for two 
isobars, namely, p = 1 bar and p = 1000 bar. The results are presented in figure [71 As can 
be seen, the model hints at the existence of minima in the heat capacity for both of these 
isobars. The location of the minimum seems to move to lower temperatures as the pressure 
increases, in concordance with experiment (see fig. 7 of Ref.— ). 

In figurelHlthe values of Cp from the simulations are compared to the experimental^ values 
for water and D2O. As one can see, the model overestimates the experimental values of the 
heat capacity of water at constant pressure. The large difference in the experimental values 
of the heat capacity of II2O and D2O demonstrate the importance of nuclear quantum effects 
in the description of the heat capacity of water. Not surprisingly, the values of TIP4P/2005 
which were obtained through classical simulations lie closer to the experimental values of 
D2O, probably reflecting the somewhat more classical behaviour of D2O with respect to 
H2O. 

Taking into account the success of TIP4P/2005 in describing most of the properties of 
water, the failure in the quantitative description of Cp (along with the important differences 
in the experimental values of D2O and H2O) points out the necessity of incorporating nuclear 
quantum effects for a quantitative description of this property. Although one could incor- 
porate quantum corrections empirically^^"^, quantum simulations (for example, via path 
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integrals) are probably the only way to correctly describe the heat capacity of water— i^. 

IV. CONCLUSIONS 

In this work molecular dynamics simulations have been undertaken for liquid water along 
several isobars. The TIP4P/2005 model has been used to describe the interaction between 
water molecules. The results obtained in relatively long runs (i.e., 40 ns) indicate that a 
maximum in density is found for the isobars considered in this work (i.e., p=l, 400, 1000 
and 1500 bar). The TMD decreases by about 1 K for each 45 bar of applied pressure. 
Thus, the TIP4P/2005 model, designed to reproduce the TMD at atmospheric pressure, 
is also able to predict the dependence of the TMD with pressure. Motivated by recent 
experimental work the behaviour of the density at low temperatures along the atmospheric 
pressure isobar was also studied, resulting in a density minimum at temperatures around 
200 K. The location and density at the minimum of the model are in very nice agreement 
with recent experimental work on the equation state of water at very low temperatures 
(obtained by confining water in narrow pores in order to prevent the nucleation of ice). Also 
in agreement with experimental measurements we have found a significant increase of the 
diffusion coefficient of super-cooled water with pressure. 

Finally, we have also computed the isothermal compressibility along the atmospheric 
pressure isobar for several water models and, once again, the results for TIP4P/2005 are 
those closest to experimental values. A minimum in the isothermal compressibility was 
found for temperatures around 37°C, once again in good agreement with the location of 
the minimum found in experiments (46.5°C) . Not only is the location of the minimum well 
described by the model, but the value at the minimum as well. Finally, we have calculated 
kt at several temperatures at a higher pressure (1000 bar) observing that TIP4P/2005 also 
exhibits a compressibility minimum. In accordance with experiment, the minimum is shifted 
towards slightly higher temperature^ when the pressure increases from 1 bar to 1000 bar. 

The results of this work indicate that the TIP4P/2005 is able to reproduce almost quan- 
titatively many of the anomalous properties of water occurring at low temperatures with 
the exception of the heat capacity. Therefore, the model can be used with confidence to ob- 
tain a better understanding of the behaviour of water in the super-cooled regime where the 
experimental determination of properties is a difficult task. It seems that for super-cooled 
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water the performance of the TIP4P/2005 model of water supersedes the performance of 
other more traditional models. 
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TABLE I: Molecular dynamics results for the TIP4P/2005 model of water along the p = 1 bar 
isobar. Only the residual part of the internal energy is given. 
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TABLE II: Molecular dynamics results for the TIP4P/2005 model of water along the p=400 bar 
isobar. Only the residual part of the internal energy is given. 



T/K p/bar p/{g/cm^) U/(kJ/mol) 
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TABLE III: Molecular dynamics results for the TIP4P/2005 model of water along the p=1000 bar 
isobar. Only the residual part of the internal energy is given. 



T/K p/bar p/{g/cm^) U/(kJ/mol) 
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TABLE IV: Molecular dynamics results for the TIP4P/2005 model of water along the p=1500 bar 
isobar. Only the residual part of the internal energy is given. 



T/K p/bar ^/(g/cm^) U/(kJ/mol) 
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TABLE V: Temperature of maximum density at different pressures. The TMD has been obtained 
by fitting the data in the proximity of the maximum to a quadratic or a cubic polynomial. The 
temperature ranges used in the fit were 260-300 K for p=l bar, 242-295 K for p=400 bar, 233-283 K 
for p=1000 bar, and 207-272 K for p=1500 bar. 
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TABLE VI: Different contributions to tfie residual internal energy of water along the T = 224.6K 
isotherm. Results were obtained from molecular dynamics simulations of the TIP4P/2005 model 
of water. The Lennard-Jones contribution C/Lennard- Jones j the Coulombic contribution C/couiombic 
and the total residual energy U are reported. All energies are given in kJ/mol. 
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TABLE VIL Isothermal compressibility kt for p = 1 bar as obtained in this work for the TIP4P, 
SPC/E and TIP5P models of water. Reported values correspond to {kt/ bar~^) xlO^. 



T/K TIP4P SPC/E TIP5P Experiment 
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TABLE VIII: Isothermal compressibility for p = 1 bar and p = 1000 bar as obtained in this 
work for the TIP4P/2005 model of water. Reported values correspond to (kt/ bar~^) xlO^. 
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FIG. 1: Density of liquid water at p = 1 bar (upper panel) and p=400 bar (lower panel) as 
obtained from Monte Carlo and molecular dynamics simulations with the TIP4P/2005 model. The 



experimental results were taken from Ref. 



54. 
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FIG. 2: Dependence of the TMD with pressure for the TIP4P/2005 model (open circles). Experi- 



mental data for the TMD (solid line) taken from Fig. 8 of Ref. 
curve of ice for the TIP4P/2005 model 
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. The dashed line is the melting 
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FIG. 3: Density of liquid water at atmospheric pressure as obtained from molecular dynamics 
simulations with the TIP4P/2005 model. For comparison, experimental data of water confined in 
narrow pores are also given^. 
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FIG. 4: Density of liquid water at atmospheric pressure as obtained from molecular dynamics sim- 
ulations with the TIP4P/2005 model (solid line and open circles). Density of ice at atmospheric 
pressure as obtained from molecular dynamics simulations with the TIP4P/2005 model (dashed 
line and open squares). 
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FIG. 5: Self-diffusion coefficient (in cm^/s) of the TIP4P 
compared to the experimental results (taken from Ref. 



/2005 model at p=l bar and p=1000 bar 
53). 
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FIG. 6: Isothermal compressibility calculated from molecular dynamics simulations. Left: Results 
for SPC/E, TIP4P and TIP5P at p=l bar using the fluctuation formula (Eq. [2]). Right: Results 
for TIP4P/2005 using the fluctuation equation (squares) and the derivative route (Eq. [H circles); 
upper results are for 1 bar and lower curves are for 1000 bar. Notice that we have deliberately used 
the same scale in both panels. For comparison, experimental dataM (full lines) are also plotted. 
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FIG. 7: Heat capacity of water at constant pressure (Cp) as obtained from simulation results of 
the TIP4P/2005 for p = 1 bar (solid line/circles) and for p = 1000 bar (dashed line/squares). The 
symbols represent the results of a first degree finite difference calculation for equation (3), and the 
curves were obtained from the results of a third degree polynomial fit (p=1000 bar ) and fourth 
degree (p=l bar) for the enthalpy. 
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FIG. 8: Heat capacity of water at constant pressure ( Cp ) as obtained from simulation results of the 
TIP4P/2005 for p = 1 bar . Experimental results for D2O and for H2O taken from Angell,Oguni 
and SichinaSS are also presented. 
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